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ABSTRACT 

Although there are a wealth of column density tracers for both the molecular and diffuse interstellar 
medium, there are very few observational studies investigating the relationship between the density 
variance (a 2 ) and the sonic Mach number (A4 S ). This is in part due to the fact that the a 2 -M. s 
relationship is derived, via MHD simulations, for the 3D density variance only, which is not a direct 
observable. We investigate the utility of a 2D column density ^\j^ Q -M-s relationship using solenoidally 

driven isothermal MHD simulations and find that the best fit follows closely the form of the 3D density 



' pipv'M-s trend but includes a scaling parameter A such that: cf n (£/£ ) = A x ln(l + b 2 A4 2 ), where 
A = 0.11 and b = 1/3. This relation is consistent with the observational data reported for the Taurus 
and IC 5146 molecular clouds with b = 0.5 and A — 0.16 and b = 0.5 and A = 0.12, respectively. 
These results open up the possibility of using the 2D column density values of a 2 for investigations of 
the relation between the sonic Mach number and the PDF variance in addition to existing PDF sonic 
Mach number relations. 

Subject headings: ISM: structure — turbulence 



1. INTRODUCTION 

Magnetohydrodynamic (MHD) turbulence is a criti- 
cal ingredient to include when considering the physics 
of the interstellar medium (ISM) . The turbulent nature 
of ISM gasses is evident from a variety of observations 
including electron density fluctuations (see Armstrong 
et al. 1994; Chepurnov & Lazarian 2010), non-thermal 
broadening of emission and absorption lines, (e.g. CO, 
HI, see Spitzer 1979; Stutzki & Guesten 1990; Heyer & 
Brunt 2004) and fractal and hierarchical structures in the 
diffuse and molecular ISM (see Elmegreen & Elmegreen 
1983; Vazequez-Semandeni 1994; Stutzki et al. 1998; 
Sanchez et al. 2005; Roman-Duval et al. 2010). A 
number of new techniques, including those studying the 
turbulence velocity spectrum (see Lazarian 2009 for a 
review) and the sonic Mach number Ai s = (l v l/C s ) 
and Alfven number M.a = (|v|/Va) of turbulenccQ (see 
Kowal et al. 2007; Burkhart et al. 2009, 2010; Burkhart, 
Lazarian, & Gaensler 2012; Esquivel & Lazarian 2010; 
Tofflemire, Burkhart & Lazarian 2011) have been devel- 
oped recently. 

An additional signature of a turbulent ISM is that 
the density and column density probability distribution 
functions (PDF) are expected to take a log-normal form 
(Vazquez-Semadeni 1994)0 Log-normal PDFs have been 
observed in multiple phases in the ISM including in 
molecular clouds (Brunt 2010), in dust extinction maps 
(Kainulainen et al. 2011; Kainulainen & Tan 2012), 
and in the diffuse ISM (Hill et al. 2008; Berkhuijsen 
& Fletcher 2008). Furthermore, the PDF was shown to 
be important for analytic models of star formation rates 
and initial mass functions (Krumholz & Mckee 2005; 

1 The sound and the Alfven speed are denoted by C B and Va, 
respectively. 

2 The relation between turbulence and the log-normal distri- 
bution can be understood as a consequence of the multiplicative 
central limit theorem assuming that individual density perturba- 
tions are independent and random. 



Hennebelle & Chabrier 2008,2011; Padoan & Nordlund 
2011). 

The PDF of both the logarithmic and linear distribu- 
tion of the gas is an important method for determin- 
ing M s . In general, the most common PDF-to-7W s pre- 
scription is to utilize numerical simulations to formulate 
an empirical relationship between the PDF moments of 
the density or column density and relate these back to 
the sonic Mach number. For example, several authors 
have suggested the turbulent sonic Mach number can be 
estimated from the calculation of the density variance 
(Padoan et al. 1997; Passot & Vazquez-Semadeni 1998; 
Beetz et al. 2008; Price, Federrath, & Brunt 2011, hence- 
forth known as PFB11) and the column density skew- 
ness/kurtosis (Kowal et al. 2007; Burkhart et al. 2009, 
2010). Other studies investigate the utility of finding M. s 
by fitting a function to the PDF and using the result- 
ing fit parameters as descriptors of turbulence (e.g. the 
Tsallis function, Esquivel & Lazarian 2010; Tofflemire, 
Burkhart, & Lazarian 2011). 

In particular, the relationship between A4 S and the 
variance of the logarithm of the density distribution as 
seen in numerical models (Padoan et al. 1997; Passot 
& Vazquez-Semadeni 1998; PFB11) generally takes the 
form: 



p/po 



b 2 Mi 



(1) 



where po is the mean value of the 3D density field, b is a 
constant of order unity and a is the standard deviation of 
the density field normalized by its mean value (i.e. p/po). 
When taking the logarithm of the normalized density 
field this relationship becomes: 



= ln(l + b 2 M 2 s ) 



(2) 



where s = ln(p/po) and a s is the standard deviation of 
the logarithm of density (not to be confused with cr p / pQ ). 
The relationships above, including the values for 6, have 
been empirically derived from MHD and hydrodynamic 
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TABLE 1 

Description of the MHD simulations 



Model A4 B ,rms Ma Plasma Resolution 



1 


8.8 


1.4 


0.05 


512 3 


2 


7.5 


0.5 


0.01 


512 3 


3 


7.3 


1.5 


0.08 


512 3 


4 


6.1 


0.5 


0.01 


512 3 


5 


5.8 


1.7 


0.17 


512 3 


6 


5.4 


0.5 


0.02 


512 3 


7 


3.6 


1.5 


0.34 


512 3 


8 


■3.7 


0.5 


0.04 


512 3 


9 


2.8 


1.7 


0.7 


512 3 


10 


2.7 


0.6 


0.1 


512 3 


11 


2.1 


1.9 


2.4 


512 3 


12 


2.2 


0.7 


0.2 


512 3 


13 


0.8 


1.7 


9.0 


512 3 


14 


0.8 


0.7 


1.5 


512 3 


15 


0.6 


1.7 


16.1 


512 3 


16 


0.6 


0.7 


2.7 


512 3 


17 


0. 1 


1.7 


36.1 


512 3 


18 


0. 1 


0.7 


6.1 


512 3 


19 


0.7 


2.7 


29.7 


256 3 


20 


2.2 


3.2 


1.1 


256 3 


21 


2.0 


2.0 


2.0 


256 3 


22 


0.7 


1.9 


14.7 


256 3 


23 


0.9 


1.9 


8.9 


256 3 


24 


6.0 


0.5 


0.01 


256 3 


25 


0.7 


0.6 


1.5 


256 3 


2(5 


2.0 


0.6 


0.2 


256 3 


27 


0.9 


0.7 


1.2 


256 3 


28 


1.0 


0.5 


0.5 


256 3 


29 


2.7 


0.4 


0.04 


256 3 


30 


1.0 


0.3 


0.2 


256 3 


31 


3.0 


0.3 


0.02 


256 3 



numerical simulations. Generally, the value of b depends 
on the driving of the turbulence in question with 6=1/3 
for solenoidal forcing and b — 1 for compressive driv- 
ing (Nordlund & Padoan 1999; Federrath, Klessen, & 
Schmidt 2008; Federrath et al. 2010). Recently, Molina 
et al. (2012) derived several variants of Equation [2] ana- 
lytically from the MHD shock jump conditions, including 
the effects of strong and weak correlations between mag- 
netic field and 3D density. 

A key limitation of these a 2 -A4 s studies is that the 
relationships are derived only for 3D density, which is not 
an observable quantity. In this letter, we investigate the 
applicability of Equations [1] and [2] for synthetic column 
density maps, in order to make these methods more easily 
applicable to observations. In this case, we define £ = 
^n(E/£o), were X is the 2D column density distribution 
available from the observations. The organization of this 
letter is as follows: In section[2]we describe our numerical 
set up. In Section [3] we discuss the a 2 -A4 s relationship 
for column density. We discuss the results in Section |4] 
followed by the conclusions in Section [5] 

2. NUMERICAL SETUP 

We generate a database of 31 3D numerical simulations 
of isothermal compressible (MHD) turbulence with reso- 
lution 512 3 and 256 3 . We use the MHD code detailed in 
Cho & Lazarian 2003 and vary the input values for the 
sonic and Alfvenic Mach number. We briefly outline the 
major points of the numerical setup. 

The code is a second-order-accurate ENO scheme 



which solves the ideal MHD equations in a periodic 
box. We drive turbulence solenoidally with energy in- 
jected on the large scales. The magnetic field consists 
of the uniform background field and a fluctuating field: 
B = B cxt + b. Initially b = 0. We stress that sim- 
ulations are scale free and all units are related to the 
turnover time and energy injection scale. We divided our 
models into two groups corresponding to sub-Alfvenic 
(B cx t _ 1-0), super- Alfvenic (B cxt < 0.1) turbulence. 
Initial values of the Alfvenic Mach number for B cxt = 0.1 
cases are « 7.0 but decrease as the magnetic field is am- 
plified due to turbulence. For each group we computed 
several models with different values of gas pressure (see 
Table [T]) falling into regimes of subsonic and supersonic. 
We also include the plasma (3, i.e. (3 = Pthermal / Pmag , 
for each of our simulations in Table [1] We create the syn- 
thetic column density maps by integrating the 3D density 
cube along a given line-of-sight (LOS) parallel or perpen- 
dicular to the mean field. We calculate the average PDFs 
for three different sight- lines in order to take into account 
LOS effects. 

3. THE COLUMN DENSITY VARIANCE- MACH NUMBER 
RELATION 

We calculate the PDFs of the synthetic column den- 
sity maps of all models listed in Table Q] for both the 
linear distribution (i.e. S/Eo) & n d taking the natural 
logarithm of the maps (i.e. £ = Zn(£/£o))- Examples 
of the PDFs of the logarithmic column density (£) can 
be found in Figure [TJ The left panel shows sub-Alfvenic 
models while the right panel shows super- Alfvenic mod- 
els. All PDFs shown here are taken from simulations 
that have resolution 512 3 . All of our synthetic column 
densities show log-normal PDFs. Visually, it is clear that 
the higher the sonic Mach number, the larger the PDF 
width. Thus, we expect the same trend found in the case 
of 3D density (i.e. variance increases with M. s ) to hold 
for 2D column density distributions. 

While it maybe the case that the variance will increase 
with M. s for column density, Equation[2]is derived for 3D 
density and will not fit the 2D distribution. We note that 
Equation [2] is not the only one of its kind found in the 
literature. For example, Lemaster & Stone (2008) used 
a three parameter fit for the density mean-A^s relation 
while Molina et al. (2012) derived a variant of Equation[2] 
from the shock jump conditions, including an additional 
parameter for plasma /?. These derivations were again 
only done for the case of 3D density fields, and not for 
observable column density maps. 

We can expect the 2D column density to have similar 
behavior to the 3D density field (i.e. will be log-normal) 
in the limit that the integration column be smaller than, 
or comparable to the correlation length of the turbulence 
(Vazquez-Semandeni & Garcia 2001). We therefore con- 
sider a form similar to Equation [21 but include a new 
scaling parameter A along with the parameter b. Thus, 
the relationship for ( = Zn(£/£o) i s: 

a 2 = A x ln{l + b 2 M 2 s ) (3) 

For a log-normal distribution, the linear and logarith- 
mic variances are related by: = exp(a 2 ) — 1, Thus 
the corresponding relation for the linear variance based 
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on Equation [3] is: 

4/E = (b 2 M 2 s + l) A -l 



(4) 
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Fig. 1— PDFs of Zn(£/£ ). The left panel shows sub-Alfvenic 
models and the right panel shows super- Alfvenic models. All mod- 
els plotted here have resolution 512^. All models show log-normal 
PDFs over a range of sonic and Alfvenic Mach number. 

The density variance- M s relationship of Equation [2] 
has only one fit parameter, and here we have chosen a 
two-parameter fit with b, and A. However, due to the 
similarity between column density and density statistics, 
we can take the value of b to be the value expected for 
solcnoidal driving, i.e. b = 1/3. In this case, we find 
the best fit A to be A = 0.11. We plot the column den- 
sity variance- M s relationship in Figure [5] for the loga- 
rithmic distribution (£, top) and the linear distribution 
(bottom) . Blue symbols are for sub-Alfvenic simulations 
and black for super- Alfvenic cases. Diamond symbols 
represent simulations with resolution for 256 3 while as- 
terisk symbols are for cases with resolution 512 3 . Error 
bars are created by taking the standard deviation of val- 
ues between different sight lines. The black dotted line 
represents Equation [3] and 2] with A — 0.11 and 6=1/3 
for the top and bottom plots, respectively. 

The column density variance- M s relationship is nar- 
rower for low values of sonic Mach number and slightly 
wider at higher sonic Mach numbers. This is due to ob- 
serving at filamentary supersonic turbulence along differ- 
ent sight-lines. Examination of Figure [5] shows that the 
column density variance- M s has a slight dependency on 
the magnetic field within the error bars, however it is not 
a overwhelming trend. This was also reported in other 
studies of PDF moments (Kowal, Lazarian, & Beresnyak 
2007). There is a general increase of the variance with 
sub-Alfvenic turbulence, which can also be visibly seen in 
Figure [TJ This is probably related to the effects studied in 
Molina et al. (2012), who found that the a^, -M. a rela- 
tionship has a dependency on the plasma /3 in the super- 
Alfvenic MHD regime and that this relation breaks down 
in the sub-Alfvenic regime due to increased anisotropy. 
In the case of column density, an additional parameter 
that will affect the observed anisotropy is the line-of-sight 
chosen. Both the level of the observed anisotropy and the 
correlations between p and B change as the line-of-sight 
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changes (Burkhart et al. 2009). In this case, it is beyond 
the scope of this letter to examine how the observed line- 
of-sight changes the a^i^-Ms relationship of Molina et 
al. (2012) when translated to 2D column density. 

We find that there is no substantial difference between 
models with 256 3 and 512 3 resolution at the currently 
studied sonic Mach numbers. This is not surprising in the 
context of other works that studied the density variance. 
Price & Federrath (2010) found that the PDFs for both 
grid and SPH codes at 256 3 and 512 3 converge. However, 
PFB11 found that, at-M s > 5 the linear variance was 
affected by numerical resolution, while logarithmic vari- 
ance is independent of numerical resolution. At A1 s > 5 
our simulations become more spread about the expected 
trend, however it is not clear if this is due to numeri- 
cal issues or line-of-sight effects. Future column density 
studies extending to higher values of M s are needed to 
confirm this. 

3.1. Observations 

With the wealth of integrated intensity and dust ex- 
tinction maps, why are there so few observational stud- 
ies investigating the relationship between a 2 -A4 s ? The 
answer, in part, is because it is not straightforward to 
calculate a p , from irj^ , which is required for the ap- 
plication of Equation [2 Brunt (2010) and PFB11 in- 
vestigated the relationship of " p i po --M-s for Taurus and 
calculated the column density variance using a combi- 
nation of 13 CO and dust extinction maps. They found 
M s ~ 17.6 ± 1.8 for Taurus using velocity dispersions. 
However, this value was calculated incorrectly and an er- 
ratum is in preparation. The corrected value of M s for 
Taurus is ps 10 (private communication with C. Brunt, 
also see Kainulainen & Tan 2012). Brunt (2010) finds 
°ii/£ Av ~ 0-7 i 0-04 (with the error coming from the 
noise variance) for the Taurus dust extinction maps taken 
from Froebrich et al. (2007). Although Brunt (2010) 
used the 13 CO and the dust extinction maps to obtain 
the underlying column density variance, we feel that the 
use of the dust maps is more reliable due to a number of 
effects (i.e. 13 CO is more susceptible to excitation effects, 
optical depth effects, and has a small dynamic range of 
density), which we detail more in Section |4j In addition 
to the Taurus molecular cloud, Padoan et al. 1997 inves- 
tigated the a 2 -M s in the cloud IC 5146 using dust ex- 
tinction observations by Lada et al. (1994). They found 
M s « 10 and cr|/ So Av = 0.49 ± 0.01, which they then 
used to convert to a 3D density variance. The values of 
b are conservatively constrained from 0.3 — 0.8 for these 
molecular clouds (Brunt 2010) with b — 0.48 — 0.5 being 
the most agreed upon value (Padoan et al. 1997; Brunt 
2010). For Taurus, the value of b is still w 0.5 despite 
the change in the sonic Mach number (see Kainulainen 
& Tan 2012). 

Because we are interested in investigating how the 
a T,/T. ~^ s relationship can be directly applied to the 
observations, we must consider the effects of instrument 
noise and telescope smoothing on our synthetic column 
density maps. We apply Gaussian white noise with 
mean signal-to-noise=100 and smoothing with a Gaus- 
sian beam to the models listed in Table [T] We choose a 
beam size of 4.6 arcminutes, with an assumed distance 
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Fig. 2. — The column density variance- M, 
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relationship (<r|, vs _M S ). Blue symbols are for sub-Alfvenic simulations and black for 



super-Alfvenic cases. Diamond symbols represent simulations with resolution for 256 3 while astrix symbols are for cases with resolution 
512 3 . Error bars are created by taking the standard deviation of values between different sight lines. The black dotted lines represent 
Equation [3] and |4] with A = 0.11 and b = 1/3 for the top and bottom plots, respectively. 



to our 'cloud' of 140pc and a cloud size of 24pc. 

We plot the column density variance derived from the 
dust extinction maps for Taurus and IC 5146 and our 
models from Table [TJ which now include 4.6 arcminutc 
smoothing and Gaussian noise, in Figure [31 in order to 
test our a 2 -A4 s relationship given in Equation [3J and [4] 
on the observations. The linear column density variance- 
's relationship is shown in the bottom panel and the 
logarithmic distribution (£) in the top panel. Blue sym- 
bols are for sub-Alfvenic simulations and black for super- 
Alfvenic simulations. Diamond symbols represent simu- 
lations with resolution for 256 3 while asterisk symbols 
are for cases with resolution 512 3 . The black lines rep- 
resent Equation [3] and 2] (for the top and bottom plots, 
respectively) with 6=1/3 and A — 0.11 for the isother- 
mal simulations and b = 0.5 and A = 0.16 and b = 0.5 
and A = 0.12 for Taurus and IC 5146, respectively. The 
green and red triangles are the 0£/£ o --A4 s values from 
the dust extinction maps for the IC 5146 and Taurus 
molecular clouds, as given by Padoan et al. (1997) and 
Brunt (2010), respectively. 

Comparison of Figure [5] (where no smoothing is in- 
cluded) with Figure [3] shows that the inclusion of smooth- 
ing and noise, at least at the level added here, does not 
severely affect the applicability of Equations [3J and H) In 
general, smoothing decreases the value of er 2 , causing one 
to underestimate Ai s . 

The values for the Taurus and IC 5146 clouds lie very 
close to what is expected from Equations |3j and [4] for 
b = 0.5 and A = 0.16 and b = 0.5 and A = 0.12, respec- 
tively. Although Equations [3] and 0] fit well, it is clear 
that the solenoidal simulations do not fit well with the 



molecular cloud data. This was also the case in PFB11, 
who attributed the discrepancy between their simula- 
tions and observations to the lack of self-gravity and/or 
use of purely solenoidal forcing. For PFB11 the assump- 
tions involved to convert the 2D to 3D variance and the 
use of CO in the variance calculation may also effect the 
relation. The fact that the dust extinction observations 
are able to fit Equations [3J and @] using the values of b 
obtained from previous studies opens up the prospects 
of using the 2D column density variance in addition to 
the assumed 3D variance. Additionally, the parameter 
A may have dependencies on other physics, which are 
discussed below. 

4. DISCUSSION 

In this letter we found a relation between the variance 
of column density and A4 S . The 3D density variance re- 
lationship of Equation [5] was modified to obtain the a 2 - 
A4 S relationship for column density. This was motivated 
by the fact that a 1 , is not available from observations. 

J PI Pa 

Ideally, both Cp/ po , as determined by the method out- 
lined in Brunt, Federrath, Price (2010), and the direct 
relationship outlined here should be used to- 
gether to provide more confidence in the a 2 -A4 s relation 
and help constrain the values of b. 

Furthermore, the results of this letter can be used syn- 
ergistically with higher order PDF moments suggested 
in a number of recent papers (see Kowal et al. 2007, 

3 Molecular clouds most likely have a mixed solenoidal and com- 
pressive forcing environment, which should yield higher values of 
b, see Federrath et al. (2010) 
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Fig. 3. — Similar to Figure[2]only now the synthetic column density maps include smoothing and Gaussian noise. The trend is unaffected 
at lower Mach number regimes, however at high sonic Mach number (M s > 7.0) the values for a 2 show slightly more of a decrease in 
variance due to shocks being smoothed over. The Taurus (red triangle; Brunt 2010) and IC 5146 (green triangle; Padoan et al. 1997) 
clouds are included along with the fits from Equations [3] and [4] 



Burkhart et al. 2009) to study M s . The use of differ- 
ent approaches in determining Ai s from observed column 
densities including other techniques, e.g. utilizing Tsal- 
lis statistics (Esquivel & Lazarian 2010, Tofflemire et al. 
2011) and the genus measure (see Burkhart, Lazarian & 
Gaensler 2012), can further increase the reliability of the 
results obtained. 

An additional complicating factor for the observations 
is that obtaining the true column density PDF can be 
tricky (Goodman, Pineda & Schnee 2009). In general, 
the use of dust extinction maps is more reliable than 
13 CO in order to obtain the true column density for a 
variety of reasons, including a larger dynamic range (fac- 
tors of 50 or more, see Goodman, Pineda & Schnee 2009, 
for more detailed discussion). Pineda, Caselli,& Good- 
man (2008) showed that the 13 CO derived column den- 
sity estimates are adversely effected by optical depth as 
low as A„ w 4 mag. Molecular transitions have a limited 
range of volume densities that they can trace, since be- 
low a critical density there may not be enough molecules 
to excite the transition and at very high densities opti- 
cal depth effects will mask the true column density. For 
these reasons, we chose to plot the values of a 2 for Tau- 
rus and IC 5146 taken from the dust extinction maps, 
rather than the 13 CO maps. 

Parameters A and b have dependencies on external 
physics beyond the sonic Mach number, b depends on 
the driving of the turbulence and it is possible that A 
may follow the same behavior. In order to test this, 
we take values from Table 3 of Federrath et al. 2010, 

with 4/S , so /=°- 46 and CT l/S ,c O m =1 - 51 and M s =5.5. 

We note that these values are taken from simulations 



that have purely solenoidal or purely compressive forc- 
ing, when in realty molecular clouds most likely have 
turbulence driven with mixed forcing. We find that ap- 
plying these values to Equation |3] yields A so i = 0.14 and 
A C om — 0.6. However, when using a simple scaling be- 
tween tJs/i; and CTp/p (given in Table 1 of Federrath et 
al. 2010), we find that A sot = 0.12 and A com = 0.25. 
Both methods show that the parameter A has a signif- 
icant dependency on the type of driving, however there 
is large variation in the compressive driven values of A 
based on these tables, which may be due to LOS effects or 
insufficient statistics. Additionally, some dependency on 
observational effects such as beam smoothing and noise 
was seen in Figure [3] and future studies will determine 
what the effects of optical depth my have on these pa- 
rameters. 

The fact of the matter is that the variance of the den- 
sity or column density distribution does not only depend 
on the sonic Mach number, but a host of other contribut- 
ing factors. This is strong motivation for the use of mul- 
tiple tools and techniques to obtain information on the 
turbulent state of the gas. A collaborative use of PDF 
methods that include the 3D and 2D variance as well as 
higher order moments of the linear distribution will yield 
the most accurate picture. 

5. CONCLUSIONS 

We investigate the utility of a 2D column density 
ct|, Eq -/V1 s trend which can be used in addition to the 



existing 3D density a 2 p , p -M. s 



trend. We find that: 



• Equations [3] and [4] (for logarithmic and linear vari- 
ance, respectively) empirically fit synthetic column 
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density maps well, even when including observa- 
tional effects such as smoothing and noise, with 
b = 1/3 (as appropriate for the solenoidal forcing 
used in this study, see Federrath et al. 2008) and 
A = 0.11. 

For the dust extinction maps of the Taurus and 
IC 5146 molecular clouds, we find A = 0.16 and 
A = 0.12, respectively, when using b = 0.5 (as 
given in the literature) for the c^^-Ads relation. 
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